{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true,
    "pycharm": {
     "name": "#%% md\n"
    }
   },
   "source": [
    "## 导入包\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "outputs": [],
   "source": [
    "# 基础的数据处理工具\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "\n",
    "# 可视化\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# 处理python时间函数\n",
    "import datetime\n",
    "\n",
    "# 处理nc数据\n",
    "import netCDF4 as nc\n",
    "from netCDF4 import num2date\n",
    "\n",
    "# 处理网格数据，shp之类的\n",
    "import geopandas as gpd\n",
    "\n",
    "# 处理tiff文件\n",
    "import rasterio\n",
    "\n",
    "# gis的一些逻辑判断\n",
    "from shapely.geometry import Point\n",
    "\n",
    "# 设置投影坐标系等\n",
    "from cartopy import crs as ccrs\n",
    "\n",
    "# 打印进度条\n",
    "from tqdm import tqdm\n",
    "\n",
    "tqdm.pandas()\n",
    "\n",
    "# 并行\n",
    "\n",
    "from joblib import Parallel, delayed\n",
    "\n",
    "# 检测系统\n",
    "\n",
    "import platform\n",
    "\n",
    "# matplotlib 显示中文的问题\n",
    "if platform.system() == 'Darwin':\n",
    "    plt.rcParams[\"font.family\"] = 'Arial Unicode MS'\n",
    "elif platform.system() == 'Windows':\n",
    "    plt.rcParams[\"font.family\"] = 'SimHei'\n",
    "else:\n",
    "    pass"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "markdown",
   "source": [
    "## 加载数据"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%% md\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "outputs": [],
   "source": [
    "shp_data = gpd.read_file(\"./数据集/Pearl王川/shp-数据1/ca_Union.shp\")\n",
    "\n",
    "nc_1988tp = nc.Dataset(\"./数据集/Pearl王川/1988tp.nc\")"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "**********************************************************************\n",
      "<class 'netCDF4._netCDF4.Variable'>\n",
      "int64 time(time)\n",
      "    units: days since 1988-01-01 00:00:00\n",
      "    calendar: proleptic_gregorian\n",
      "unlimited dimensions: \n",
      "current shape = (366,)\n",
      "filling on, default _FillValue of -9223372036854775806 used\n",
      "**********************************************************************\n",
      "<class 'netCDF4._netCDF4.Variable'>\n",
      "float32 longitude(longitude)\n",
      "    _FillValue: nan\n",
      "    units: degrees_east\n",
      "    long_name: longitude\n",
      "unlimited dimensions: \n",
      "current shape = (611,)\n",
      "filling on\n",
      "**********************************************************************\n",
      "<class 'netCDF4._netCDF4.Variable'>\n",
      "float32 latitude(latitude)\n",
      "    _FillValue: nan\n",
      "    units: degrees_north\n",
      "    long_name: latitude\n",
      "unlimited dimensions: \n",
      "current shape = (221,)\n",
      "filling on\n",
      "**********************************************************************\n",
      "<class 'netCDF4._netCDF4.Variable'>\n",
      "float32 tp(time, latitude, longitude)\n",
      "    _FillValue: nan\n",
      "unlimited dimensions: \n",
      "current shape = (366, 221, 611)\n",
      "filling on\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "<ipython-input-3-e3e8d12b5fc0>:3: DeprecationWarning: tostring() is deprecated. Use tobytes() instead.\n",
      "  print(item)\n"
     ]
    }
   ],
   "source": [
    "for item in nc_1988tp.variables.values():\n",
    "    print('*' * 70)\n",
    "    print(item)"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "<ipython-input-4-3db7c478fb2a>:1: DeprecationWarning: `np.bool` is a deprecated alias for the builtin `bool`. To silence this warning, use `bool` by itself. Doing this will not modify any behavior and is safe. If you specifically wanted the numpy scalar type, use `np.bool_` here.\n",
      "Deprecated in NumPy 1.20; for more details and guidance: https://numpy.org/devdocs/release/1.20.0-notes.html#deprecations\n",
      "  raw_longitude = np.array(nc_1988tp.variables['longitude'])\n",
      "<ipython-input-4-3db7c478fb2a>:2: DeprecationWarning: `np.bool` is a deprecated alias for the builtin `bool`. To silence this warning, use `bool` by itself. Doing this will not modify any behavior and is safe. If you specifically wanted the numpy scalar type, use `np.bool_` here.\n",
      "Deprecated in NumPy 1.20; for more details and guidance: https://numpy.org/devdocs/release/1.20.0-notes.html#deprecations\n",
      "  raw_latitude = np.array(nc_1988tp.variables['latitude'])\n",
      "<ipython-input-4-3db7c478fb2a>:3: DeprecationWarning: `np.bool` is a deprecated alias for the builtin `bool`. To silence this warning, use `bool` by itself. Doing this will not modify any behavior and is safe. If you specifically wanted the numpy scalar type, use `np.bool_` here.\n",
      "Deprecated in NumPy 1.20; for more details and guidance: https://numpy.org/devdocs/release/1.20.0-notes.html#deprecations\n",
      "  raw_time = np.array(nc_1988tp.variables['time'])\n",
      "<ipython-input-4-3db7c478fb2a>:4: DeprecationWarning: `np.bool` is a deprecated alias for the builtin `bool`. To silence this warning, use `bool` by itself. Doing this will not modify any behavior and is safe. If you specifically wanted the numpy scalar type, use `np.bool_` here.\n",
      "Deprecated in NumPy 1.20; for more details and guidance: https://numpy.org/devdocs/release/1.20.0-notes.html#deprecations\n",
      "  raw_tp = np.array(nc_1988tp.variables['tp'])\n"
     ]
    }
   ],
   "source": [
    "raw_longitude = np.array(nc_1988tp.variables['longitude'])\n",
    "raw_latitude = np.array(nc_1988tp.variables['latitude'])\n",
    "raw_time = np.array(nc_1988tp.variables['time'])\n",
    "raw_tp = np.array(nc_1988tp.variables['tp'])"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "markdown",
   "source": [
    "## 插值"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%% md\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "100%|██████████| 366/366 [00:02<00:00, 140.48it/s]\n"
     ]
    },
    {
     "data": {
      "text/plain": "(366, 360, 720)"
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# [ 89.75 -89.75]\n",
    "# [-179.75 179.75]\n",
    "target_lon = -179.75 + 0.5 * np.arange(0, 720)\n",
    "target_lat = -89.75 + 0.5 * np.arange(0, 360)\n",
    "\n",
    "target_value = []\n",
    "\n",
    "from scipy.interpolate import interp2d\n",
    "\n",
    "for i in tqdm(range(raw_time.shape[0])):\n",
    "    f = interp2d(raw_longitude, raw_latitude, raw_tp[i, :, :])\n",
    "    temp_value = f(target_lon, target_lat)\n",
    "    target_value.append(temp_value)\n",
    "\n",
    "target_value = np.array(target_value)\n",
    "target_value.shape"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "markdown",
   "source": [
    "## 写入nc文件"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%% md\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "outputs": [],
   "source": [
    "with nc.Dataset('./结果/test1115.nc', mode='w', format='NETCDF4_CLASSIC') as ncfile:\n",
    "    # 创建维度\n",
    "    lat_dim = ncfile.createDimension('lat', 360)  # latitude axis\n",
    "    lon_dim = ncfile.createDimension('lon', 720)  # longitude axis\n",
    "    time_dim = ncfile.createDimension('time', None)\n",
    "\n",
    "    # 创建变量\n",
    "    lat = ncfile.createVariable('lat', np.float32, ('lat',))\n",
    "    lat.units = 'degrees_north'\n",
    "    lat.long_name = 'latitude'\n",
    "\n",
    "    lon = ncfile.createVariable('lon', np.float32, ('lon',))\n",
    "    lon.units = 'degrees_east'\n",
    "    lon.long_name = 'longitude'\n",
    "\n",
    "    time = ncfile.createVariable('time', np.float64, ('time',))\n",
    "    time.units = 'days since 1988-01-01 00:00:00'\n",
    "    time.long_name = 'time'\n",
    "\n",
    "    temp = ncfile.createVariable('temp', np.float64, ('time', 'lat', 'lon'))  # note: unlimited dimension is leftmost\n",
    "    # temp.units = 'K' # degrees Kelvin\n",
    "    temp.standard_name = 'air_temperature'\n",
    "\n",
    "    # 写入变量\n",
    "    lat[:] = target_lat\n",
    "    lon[:] = target_lon\n",
    "    time[:] = raw_time\n",
    "    temp[:, :, :] = target_value"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "markdown",
   "source": [
    "### 测试"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%% md\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "outputs": [
    {
     "data": {
      "text/plain": "<xarray.Dataset>\nDimensions:  (lat: 360, lon: 720, time: 366)\nCoordinates:\n  * lat      (lat) float32 -89.75 -89.25 -88.75 -88.25 ... 88.75 89.25 89.75\n  * lon      (lon) float32 -179.8 -179.2 -178.8 -178.2 ... 178.8 179.2 179.8\n  * time     (time) datetime64[ns] 1988-01-01 1988-01-02 ... 1988-12-31\nData variables:\n    temp     (time, lat, lon) float64 ...",
      "text/html": "<div><svg style=\"position: absolute; width: 0; height: 0; overflow: hidden\">\n<defs>\n<symbol id=\"icon-database\" viewBox=\"0 0 32 32\">\n<path d=\"M16 0c-8.837 0-16 2.239-16 5v4c0 2.761 7.163 5 16 5s16-2.239 16-5v-4c0-2.761-7.163-5-16-5z\"></path>\n<path d=\"M16 17c-8.837 0-16-2.239-16-5v6c0 2.761 7.163 5 16 5s16-2.239 16-5v-6c0 2.761-7.163 5-16 5z\"></path>\n<path d=\"M16 26c-8.837 0-16-2.239-16-5v6c0 2.761 7.163 5 16 5s16-2.239 16-5v-6c0 2.761-7.163 5-16 5z\"></path>\n</symbol>\n<symbol id=\"icon-file-text2\" viewBox=\"0 0 32 32\">\n<path d=\"M28.681 7.159c-0.694-0.947-1.662-2.053-2.724-3.116s-2.169-2.030-3.116-2.724c-1.612-1.182-2.393-1.319-2.841-1.319h-15.5c-1.378 0-2.5 1.121-2.5 2.5v27c0 1.378 1.122 2.5 2.5 2.5h23c1.378 0 2.5-1.122 2.5-2.5v-19.5c0-0.448-0.137-1.23-1.319-2.841zM24.543 5.457c0.959 0.959 1.712 1.825 2.268 2.543h-4.811v-4.811c0.718 0.556 1.584 1.309 2.543 2.268zM28 29.5c0 0.271-0.229 0.5-0.5 0.5h-23c-0.271 0-0.5-0.229-0.5-0.5v-27c0-0.271 0.229-0.5 0.5-0.5 0 0 15.499-0 15.5 0v7c0 0.552 0.448 1 1 1h7v19.5z\"></path>\n<path d=\"M23 26h-14c-0.552 0-1-0.448-1-1s0.448-1 1-1h14c0.552 0 1 0.448 1 1s-0.448 1-1 1z\"></path>\n<path d=\"M23 22h-14c-0.552 0-1-0.448-1-1s0.448-1 1-1h14c0.552 0 1 0.448 1 1s-0.448 1-1 1z\"></path>\n<path d=\"M23 18h-14c-0.552 0-1-0.448-1-1s0.448-1 1-1h14c0.552 0 1 0.448 1 1s-0.448 1-1 1z\"></path>\n</symbol>\n</defs>\n</svg>\n<style>/* CSS stylesheet for displaying xarray objects in jupyterlab.\n *\n */\n\n:root {\n  --xr-font-color0: var(--jp-content-font-color0, rgba(0, 0, 0, 1));\n  --xr-font-color2: var(--jp-content-font-color2, rgba(0, 0, 0, 0.54));\n  --xr-font-color3: var(--jp-content-font-color3, rgba(0, 0, 0, 0.38));\n  --xr-border-color: var(--jp-border-color2, #e0e0e0);\n  --xr-disabled-color: var(--jp-layout-color3, #bdbdbd);\n  --xr-background-color: var(--jp-layout-color0, white);\n  --xr-background-color-row-even: var(--jp-layout-color1, white);\n  --xr-background-color-row-odd: var(--jp-layout-color2, #eeeeee);\n}\n\nhtml[theme=dark],\nbody.vscode-dark {\n  --xr-font-color0: rgba(255, 255, 255, 1);\n  --xr-font-color2: rgba(255, 255, 255, 0.54);\n  --xr-font-color3: rgba(255, 255, 255, 0.38);\n  --xr-border-color: #1F1F1F;\n  --xr-disabled-color: #515151;\n  --xr-background-color: #111111;\n  --xr-background-color-row-even: #111111;\n  --xr-background-color-row-odd: #313131;\n}\n\n.xr-wrap {\n  display: block;\n  min-width: 300px;\n  max-width: 700px;\n}\n\n.xr-text-repr-fallback {\n  /* fallback to plain text repr when CSS is not injected (untrusted notebook) */\n  display: none;\n}\n\n.xr-header {\n  padding-top: 6px;\n  padding-bottom: 6px;\n  margin-bottom: 4px;\n  border-bottom: solid 1px var(--xr-border-color);\n}\n\n.xr-header > div,\n.xr-header > ul {\n  display: inline;\n  margin-top: 0;\n  margin-bottom: 0;\n}\n\n.xr-obj-type,\n.xr-array-name {\n  margin-left: 2px;\n  margin-right: 10px;\n}\n\n.xr-obj-type {\n  color: var(--xr-font-color2);\n}\n\n.xr-sections {\n  padding-left: 0 !important;\n  display: grid;\n  grid-template-columns: 150px auto auto 1fr 20px 20px;\n}\n\n.xr-section-item {\n  display: contents;\n}\n\n.xr-section-item input {\n  display: none;\n}\n\n.xr-section-item input + label {\n  color: var(--xr-disabled-color);\n}\n\n.xr-section-item input:enabled + label {\n  cursor: pointer;\n  color: var(--xr-font-color2);\n}\n\n.xr-section-item input:enabled + label:hover {\n  color: var(--xr-font-color0);\n}\n\n.xr-section-summary {\n  grid-column: 1;\n  color: var(--xr-font-color2);\n  font-weight: 500;\n}\n\n.xr-section-summary > span {\n  display: inline-block;\n  padding-left: 0.5em;\n}\n\n.xr-section-summary-in:disabled + label {\n  color: var(--xr-font-color2);\n}\n\n.xr-section-summary-in + label:before {\n  display: inline-block;\n  content: '►';\n  font-size: 11px;\n  width: 15px;\n  text-align: center;\n}\n\n.xr-section-summary-in:disabled + label:before {\n  color: var(--xr-disabled-color);\n}\n\n.xr-section-summary-in:checked + label:before {\n  content: '▼';\n}\n\n.xr-section-summary-in:checked + label > span {\n  display: none;\n}\n\n.xr-section-summary,\n.xr-section-inline-details {\n  padding-top: 4px;\n  padding-bottom: 4px;\n}\n\n.xr-section-inline-details {\n  grid-column: 2 / -1;\n}\n\n.xr-section-details {\n  display: none;\n  grid-column: 1 / -1;\n  margin-bottom: 5px;\n}\n\n.xr-section-summary-in:checked ~ .xr-section-details {\n  display: contents;\n}\n\n.xr-array-wrap {\n  grid-column: 1 / -1;\n  display: grid;\n  grid-template-columns: 20px auto;\n}\n\n.xr-array-wrap > label {\n  grid-column: 1;\n  vertical-align: top;\n}\n\n.xr-preview {\n  color: var(--xr-font-color3);\n}\n\n.xr-array-preview,\n.xr-array-data {\n  padding: 0 5px !important;\n  grid-column: 2;\n}\n\n.xr-array-data,\n.xr-array-in:checked ~ .xr-array-preview {\n  display: none;\n}\n\n.xr-array-in:checked ~ .xr-array-data,\n.xr-array-preview {\n  display: inline-block;\n}\n\n.xr-dim-list {\n  display: inline-block !important;\n  list-style: none;\n  padding: 0 !important;\n  margin: 0;\n}\n\n.xr-dim-list li {\n  display: inline-block;\n  padding: 0;\n  margin: 0;\n}\n\n.xr-dim-list:before {\n  content: '(';\n}\n\n.xr-dim-list:after {\n  content: ')';\n}\n\n.xr-dim-list li:not(:last-child):after {\n  content: ',';\n  padding-right: 5px;\n}\n\n.xr-has-index {\n  font-weight: bold;\n}\n\n.xr-var-list,\n.xr-var-item {\n  display: contents;\n}\n\n.xr-var-item > div,\n.xr-var-item label,\n.xr-var-item > .xr-var-name span {\n  background-color: var(--xr-background-color-row-even);\n  margin-bottom: 0;\n}\n\n.xr-var-item > .xr-var-name:hover span {\n  padding-right: 5px;\n}\n\n.xr-var-list > li:nth-child(odd) > div,\n.xr-var-list > li:nth-child(odd) > label,\n.xr-var-list > li:nth-child(odd) > .xr-var-name span {\n  background-color: var(--xr-background-color-row-odd);\n}\n\n.xr-var-name {\n  grid-column: 1;\n}\n\n.xr-var-dims {\n  grid-column: 2;\n}\n\n.xr-var-dtype {\n  grid-column: 3;\n  text-align: right;\n  color: var(--xr-font-color2);\n}\n\n.xr-var-preview {\n  grid-column: 4;\n}\n\n.xr-var-name,\n.xr-var-dims,\n.xr-var-dtype,\n.xr-preview,\n.xr-attrs dt {\n  white-space: nowrap;\n  overflow: hidden;\n  text-overflow: ellipsis;\n  padding-right: 10px;\n}\n\n.xr-var-name:hover,\n.xr-var-dims:hover,\n.xr-var-dtype:hover,\n.xr-attrs dt:hover {\n  overflow: visible;\n  width: auto;\n  z-index: 1;\n}\n\n.xr-var-attrs,\n.xr-var-data {\n  display: none;\n  background-color: var(--xr-background-color) !important;\n  padding-bottom: 5px !important;\n}\n\n.xr-var-attrs-in:checked ~ .xr-var-attrs,\n.xr-var-data-in:checked ~ .xr-var-data {\n  display: block;\n}\n\n.xr-var-data > table {\n  float: right;\n}\n\n.xr-var-name span,\n.xr-var-data,\n.xr-attrs {\n  padding-left: 25px !important;\n}\n\n.xr-attrs,\n.xr-var-attrs,\n.xr-var-data {\n  grid-column: 1 / -1;\n}\n\ndl.xr-attrs {\n  padding: 0;\n  margin: 0;\n  display: grid;\n  grid-template-columns: 125px auto;\n}\n\n.xr-attrs dt,\n.xr-attrs dd {\n  padding: 0;\n  margin: 0;\n  float: left;\n  padding-right: 10px;\n  width: auto;\n}\n\n.xr-attrs dt {\n  font-weight: normal;\n  grid-column: 1;\n}\n\n.xr-attrs dt:hover span {\n  display: inline-block;\n  background: var(--xr-background-color);\n  padding-right: 10px;\n}\n\n.xr-attrs dd {\n  grid-column: 2;\n  white-space: pre-wrap;\n  word-break: break-all;\n}\n\n.xr-icon-database,\n.xr-icon-file-text2 {\n  display: inline-block;\n  vertical-align: middle;\n  width: 1em;\n  height: 1.5em !important;\n  stroke-width: 0;\n  stroke: currentColor;\n  fill: currentColor;\n}\n</style><pre class='xr-text-repr-fallback'>&lt;xarray.Dataset&gt;\nDimensions:  (lat: 360, lon: 720, time: 366)\nCoordinates:\n  * lat      (lat) float32 -89.75 -89.25 -88.75 -88.25 ... 88.75 89.25 89.75\n  * lon      (lon) float32 -179.8 -179.2 -178.8 -178.2 ... 178.8 179.2 179.8\n  * time     (time) datetime64[ns] 1988-01-01 1988-01-02 ... 1988-12-31\nData variables:\n    temp     (time, lat, lon) float64 ...</pre><div class='xr-wrap' hidden><div class='xr-header'><div class='xr-obj-type'>xarray.Dataset</div></div><ul class='xr-sections'><li class='xr-section-item'><input id='section-ab9fe715-a7a8-4d36-8a10-062200c4c0ae' class='xr-section-summary-in' type='checkbox' disabled ><label for='section-ab9fe715-a7a8-4d36-8a10-062200c4c0ae' class='xr-section-summary'  title='Expand/collapse section'>Dimensions:</label><div class='xr-section-inline-details'><ul class='xr-dim-list'><li><span class='xr-has-index'>lat</span>: 360</li><li><span class='xr-has-index'>lon</span>: 720</li><li><span class='xr-has-index'>time</span>: 366</li></ul></div><div class='xr-section-details'></div></li><li class='xr-section-item'><input id='section-a215cf74-c648-4fe4-9a89-eda997acbc27' class='xr-section-summary-in' type='checkbox'  checked><label for='section-a215cf74-c648-4fe4-9a89-eda997acbc27' class='xr-section-summary' >Coordinates: <span>(3)</span></label><div class='xr-section-inline-details'></div><div class='xr-section-details'><ul class='xr-var-list'><li class='xr-var-item'><div class='xr-var-name'><span class='xr-has-index'>lat</span></div><div class='xr-var-dims'>(lat)</div><div class='xr-var-dtype'>float32</div><div class='xr-var-preview xr-preview'>-89.75 -89.25 ... 89.25 89.75</div><input id='attrs-8e241cdd-97df-4052-abff-4e8fb406df50' class='xr-var-attrs-in' type='checkbox' ><label for='attrs-8e241cdd-97df-4052-abff-4e8fb406df50' title='Show/Hide attributes'><svg class='icon xr-icon-file-text2'><use xlink:href='#icon-file-text2'></use></svg></label><input id='data-b03f724a-f8a8-4b47-92d8-91ba8e20f371' class='xr-var-data-in' type='checkbox'><label for='data-b03f724a-f8a8-4b47-92d8-91ba8e20f371' title='Show/Hide data repr'><svg class='icon xr-icon-database'><use xlink:href='#icon-database'></use></svg></label><div class='xr-var-attrs'><dl class='xr-attrs'><dt><span>units :</span></dt><dd>degrees_north</dd><dt><span>long_name :</span></dt><dd>latitude</dd></dl></div><div class='xr-var-data'><pre>array([-89.75, -89.25, -88.75, ...,  88.75,  89.25,  89.75], dtype=float32)</pre></div></li><li class='xr-var-item'><div class='xr-var-name'><span class='xr-has-index'>lon</span></div><div class='xr-var-dims'>(lon)</div><div class='xr-var-dtype'>float32</div><div class='xr-var-preview xr-preview'>-179.8 -179.2 ... 179.2 179.8</div><input id='attrs-e545de49-b4b3-40fa-a91e-a5b2fbcc29c1' class='xr-var-attrs-in' type='checkbox' ><label for='attrs-e545de49-b4b3-40fa-a91e-a5b2fbcc29c1' title='Show/Hide attributes'><svg class='icon xr-icon-file-text2'><use xlink:href='#icon-file-text2'></use></svg></label><input id='data-877b1343-1b5f-4563-b07e-27548d543c4a' class='xr-var-data-in' type='checkbox'><label for='data-877b1343-1b5f-4563-b07e-27548d543c4a' title='Show/Hide data repr'><svg class='icon xr-icon-database'><use xlink:href='#icon-database'></use></svg></label><div class='xr-var-attrs'><dl class='xr-attrs'><dt><span>units :</span></dt><dd>degrees_east</dd><dt><span>long_name :</span></dt><dd>longitude</dd></dl></div><div class='xr-var-data'><pre>array([-179.75, -179.25, -178.75, ...,  178.75,  179.25,  179.75],\n      dtype=float32)</pre></div></li><li class='xr-var-item'><div class='xr-var-name'><span class='xr-has-index'>time</span></div><div class='xr-var-dims'>(time)</div><div class='xr-var-dtype'>datetime64[ns]</div><div class='xr-var-preview xr-preview'>1988-01-01 ... 1988-12-31</div><input id='attrs-c40566c7-b73f-4199-ae95-b7b5e9cc0c1a' class='xr-var-attrs-in' type='checkbox' ><label for='attrs-c40566c7-b73f-4199-ae95-b7b5e9cc0c1a' title='Show/Hide attributes'><svg class='icon xr-icon-file-text2'><use xlink:href='#icon-file-text2'></use></svg></label><input id='data-2cf5bf8c-5f80-44d7-b6b1-fbc809f15d40' class='xr-var-data-in' type='checkbox'><label for='data-2cf5bf8c-5f80-44d7-b6b1-fbc809f15d40' title='Show/Hide data repr'><svg class='icon xr-icon-database'><use xlink:href='#icon-database'></use></svg></label><div class='xr-var-attrs'><dl class='xr-attrs'><dt><span>long_name :</span></dt><dd>time</dd></dl></div><div class='xr-var-data'><pre>array([&#x27;1988-01-01T00:00:00.000000000&#x27;, &#x27;1988-01-02T00:00:00.000000000&#x27;,\n       &#x27;1988-01-03T00:00:00.000000000&#x27;, ..., &#x27;1988-12-29T00:00:00.000000000&#x27;,\n       &#x27;1988-12-30T00:00:00.000000000&#x27;, &#x27;1988-12-31T00:00:00.000000000&#x27;],\n      dtype=&#x27;datetime64[ns]&#x27;)</pre></div></li></ul></div></li><li class='xr-section-item'><input id='section-454652a9-3b02-434f-b051-78f16ab93b3a' class='xr-section-summary-in' type='checkbox'  checked><label for='section-454652a9-3b02-434f-b051-78f16ab93b3a' class='xr-section-summary' >Data variables: <span>(1)</span></label><div class='xr-section-inline-details'></div><div class='xr-section-details'><ul class='xr-var-list'><li class='xr-var-item'><div class='xr-var-name'><span>temp</span></div><div class='xr-var-dims'>(time, lat, lon)</div><div class='xr-var-dtype'>float64</div><div class='xr-var-preview xr-preview'>...</div><input id='attrs-2f39e07e-8604-4b8e-8167-69ac2ee26be8' class='xr-var-attrs-in' type='checkbox' ><label for='attrs-2f39e07e-8604-4b8e-8167-69ac2ee26be8' title='Show/Hide attributes'><svg class='icon xr-icon-file-text2'><use xlink:href='#icon-file-text2'></use></svg></label><input id='data-ebbf7749-f89c-4853-a9a2-82202f98a6c5' class='xr-var-data-in' type='checkbox'><label for='data-ebbf7749-f89c-4853-a9a2-82202f98a6c5' title='Show/Hide data repr'><svg class='icon xr-icon-database'><use xlink:href='#icon-database'></use></svg></label><div class='xr-var-attrs'><dl class='xr-attrs'><dt><span>standard_name :</span></dt><dd>air_temperature</dd></dl></div><div class='xr-var-data'><pre>[94867200 values with dtype=float64]</pre></div></li></ul></div></li><li class='xr-section-item'><input id='section-e619281b-a1f6-44e5-9546-fde4c7fb844e' class='xr-section-summary-in' type='checkbox' disabled ><label for='section-e619281b-a1f6-44e5-9546-fde4c7fb844e' class='xr-section-summary'  title='Expand/collapse section'>Attributes: <span>(0)</span></label><div class='xr-section-inline-details'></div><div class='xr-section-details'><dl class='xr-attrs'></dl></div></li></ul></div></div>"
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import xarray as xr\n",
    "\n",
    "test_xr = xr.open_dataset(\"./结果/test1115.nc\")\n",
    "test_xr"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "**********************************************************************\n",
      "lat\n",
      "<class 'netCDF4._netCDF4.Variable'>\n",
      "float32 lat(lat)\n",
      "    units: degrees_north\n",
      "    long_name: latitude\n",
      "unlimited dimensions: \n",
      "current shape = (360,)\n",
      "filling on, default _FillValue of 9.969209968386869e+36 used\n",
      "**********************************************************************\n",
      "lon\n",
      "<class 'netCDF4._netCDF4.Variable'>\n",
      "float32 lon(lon)\n",
      "    units: degrees_east\n",
      "    long_name: longitude\n",
      "unlimited dimensions: \n",
      "current shape = (720,)\n",
      "filling on, default _FillValue of 9.969209968386869e+36 used\n",
      "**********************************************************************\n",
      "time\n",
      "<class 'netCDF4._netCDF4.Variable'>\n",
      "float64 time(time)\n",
      "    units: days since 1988-01-01 00:00:00\n",
      "    long_name: time\n",
      "unlimited dimensions: time\n",
      "current shape = (366,)\n",
      "filling on, default _FillValue of 9.969209968386869e+36 used\n",
      "**********************************************************************\n",
      "temp\n",
      "<class 'netCDF4._netCDF4.Variable'>\n",
      "float64 temp(time, lat, lon)\n",
      "    standard_name: air_temperature\n",
      "unlimited dimensions: time\n",
      "current shape = (366, 360, 720)\n",
      "filling on, default _FillValue of 9.969209968386869e+36 used\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "<ipython-input-43-32037f2bb0f0>:6: DeprecationWarning: tostring() is deprecated. Use tobytes() instead.\n",
      "  print(item)\n"
     ]
    }
   ],
   "source": [
    "test_nc = nc.Dataset(\"./结果/test1115.nc\")\n",
    "\n",
    "for item in test_nc.variables.values():\n",
    "    print('*' * 70)\n",
    "    print(item.name)\n",
    "    print(item)"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "markdown",
   "source": [
    "## 裁切"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%% md\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "<ipython-input-44-1cb7bc6e1b9d>:76: DeprecationWarning: `np.bool` is a deprecated alias for the builtin `bool`. To silence this warning, use `bool` by itself. Doing this will not modify any behavior and is safe. If you specifically wanted the numpy scalar type, use `np.bool_` here.\n",
      "Deprecated in NumPy 1.20; for more details and guidance: https://numpy.org/devdocs/release/1.20.0-notes.html#deprecations\n",
      "  self.nc_target_data = np.array(self.nc_data.variables[self.nc_variable])\n",
      "<ipython-input-44-1cb7bc6e1b9d>:78: DeprecationWarning: tostring() is deprecated. Use tobytes() instead.\n",
      "  if 'missing_value' in dir(self.nc_data.variables[self.nc_variable]):\n",
      "<ipython-input-44-1cb7bc6e1b9d>:86: DeprecationWarning: `np.bool` is a deprecated alias for the builtin `bool`. To silence this warning, use `bool` by itself. Doing this will not modify any behavior and is safe. If you specifically wanted the numpy scalar type, use `np.bool_` here.\n",
      "Deprecated in NumPy 1.20; for more details and guidance: https://numpy.org/devdocs/release/1.20.0-notes.html#deprecations\n",
      "  self.longitude_data = np.array(self.nc_data.variables[self.lon_variable])\n",
      "<ipython-input-44-1cb7bc6e1b9d>:87: DeprecationWarning: `np.bool` is a deprecated alias for the builtin `bool`. To silence this warning, use `bool` by itself. Doing this will not modify any behavior and is safe. If you specifically wanted the numpy scalar type, use `np.bool_` here.\n",
      "Deprecated in NumPy 1.20; for more details and guidance: https://numpy.org/devdocs/release/1.20.0-notes.html#deprecations\n",
      "  self.latitude_data = np.array(self.nc_data.variables[self.lat_variable])\n",
      "<ipython-input-44-1cb7bc6e1b9d>:88: DeprecationWarning: tostring() is deprecated. Use tobytes() instead.\n",
      "  self.time_units = self.nc_data.variables[self.time_variable].units\n",
      "<ipython-input-44-1cb7bc6e1b9d>:89: DeprecationWarning: `np.bool` is a deprecated alias for the builtin `bool`. To silence this warning, use `bool` by itself. Doing this will not modify any behavior and is safe. If you specifically wanted the numpy scalar type, use `np.bool_` here.\n",
      "Deprecated in NumPy 1.20; for more details and guidance: https://numpy.org/devdocs/release/1.20.0-notes.html#deprecations\n",
      "  self.time_data = np.array(self.nc_data.variables[self.time_variable])\n",
      "100%|██████████| 720/720 [00:43<00:00, 16.72it/s]\n"
     ]
    }
   ],
   "source": [
    "class GetMask(object):\n",
    "    def __init__(self,\n",
    "                 geopandas: gpd.GeoDataFrame,\n",
    "                 nc_data: nc.Dataset,\n",
    "                 nc_variable: str,\n",
    "                 lat_variable: str,\n",
    "                 lon_variable: str,\n",
    "                 time_variable: str):\n",
    "        self.geopandas = geopandas\n",
    "        self.nc_data = nc_data\n",
    "        self.nc_variable = nc_variable\n",
    "        self.lat_variable = lat_variable\n",
    "        self.lon_variable = lon_variable\n",
    "        self.time_variable = time_variable\n",
    "        self.nc_target_data = None\n",
    "        self.target_data_missing_value = None\n",
    "        self.time_dim = None\n",
    "        self.lat_dim = None\n",
    "        self.lon_dim = None\n",
    "        self.mask_matrix = None\n",
    "        self.longitude_data = None\n",
    "        self.latitude_data = None\n",
    "        self.time_data = None\n",
    "        self.time_units = None\n",
    "        self.clean_time_data = None\n",
    "\n",
    "    def num2datetime(self, cftime, units, format='%Y-%m-%d %H:%M:%S'):\n",
    "        \"\"\"\n",
    "        将nc文件里面的时间格式 从cftime 转换到 datetime格式\n",
    "        :param cftime:\n",
    "        :param units:\n",
    "        :param format:\n",
    "        :return:\n",
    "        \"\"\"\n",
    "        return datetime.datetime.strptime(num2date(times=cftime, units=units).strftime(format), format)\n",
    "\n",
    "    @staticmethod\n",
    "    def array2gtiff(array, filename):\n",
    "        \"\"\"\n",
    "        将一个矩阵保存为tiff文件,\n",
    "        这里还可以设置tiff的crs和transofrm。更多，可以查看rasterio的官网或者下面的这个链接\n",
    "        https://gis.stackexchange.com/questions/279953/numpy-array-to-gtiff-using-rasterio-without-source-raster\n",
    "        :param array: shape:(row, col)\n",
    "        :param filename:\n",
    "        :return:\n",
    "        \"\"\"\n",
    "        with rasterio.open(filename, 'w', driver='GTiff',\n",
    "                           height=array.shape[0], width=array.shape[1],\n",
    "                           count=1, dtype=str(array.dtype)) as f:\n",
    "            f.write(array, 1)\n",
    "\n",
    "    def pic(self, lon, lat) -> bool:\n",
    "\n",
    "        \"\"\"\n",
    "        检测一个点是否在中国边界线内\n",
    "        lon:东经\n",
    "        lat:北纬\n",
    "        :param lon:\n",
    "        :param lat:\n",
    "        :return:\n",
    "        \"\"\"\n",
    "        return self.geopandas.contains(Point(lon, lat))[0]\n",
    "\n",
    "    def parallel_mask(self, index_lon, index_lat):\n",
    "        point = (self.longitude_data[index_lon], self.latitude_data[index_lat])\n",
    "        value = self.pic(lon=point[0], lat=point[1])\n",
    "        # return value\n",
    "        self.mask_matrix[index_lat, index_lon] = value\n",
    "\n",
    "    def run(self):\n",
    "        # 处理geopandas数据\n",
    "        # self.geopandas = self.geopandas.iloc[0, :]\n",
    "        self.geopandas['geometry'] = self.geopandas.buffer(0)\n",
    "\n",
    "        # 处理nc数据\n",
    "        self.nc_target_data = np.array(self.nc_data.variables[self.nc_variable])\n",
    "\n",
    "        if 'missing_value' in dir(self.nc_data.variables[self.nc_variable]):\n",
    "            self.target_data_missing_value = self.nc_data.variables[self.nc_variable].missing_value\n",
    "        else:\n",
    "            self.target_data_missing_value = np.nan\n",
    "\n",
    "        self.nc_target_data[self.nc_target_data == self.target_data_missing_value] = np.nan\n",
    "\n",
    "        # 提取lat,lon,lat 变量\n",
    "        self.longitude_data = np.array(self.nc_data.variables[self.lon_variable])\n",
    "        self.latitude_data = np.array(self.nc_data.variables[self.lat_variable])\n",
    "        self.time_units = self.nc_data.variables[self.time_variable].units\n",
    "        self.time_data = np.array(self.nc_data.variables[self.time_variable])\n",
    "        self.clean_time_data = [self.num2datetime(cftime=i, units=self.time_units) for i in self.time_data]\n",
    "\n",
    "        # 创建一个mask\n",
    "        nc_target_data_shape = self.nc_target_data.shape\n",
    "\n",
    "        if len(nc_target_data_shape) == 3:\n",
    "            (self.time_dim, self.lat_dim, self.lon_dim) = nc_target_data_shape\n",
    "        else:\n",
    "            (self.lat_dim, self.lon_dim) = nc_target_data_shape\n",
    "\n",
    "        self.mask_matrix = np.full(shape=(self.lat_dim, self.lon_dim), fill_value=False)\n",
    "\n",
    "        _ = Parallel(n_jobs=-1, backend='threading', verbose=0)(\n",
    "            delayed(self.parallel_mask)(index_lon, index_lat)\n",
    "            for index_lon in tqdm(range(self.lon_dim))\n",
    "            for index_lat in range(self.lat_dim))\n",
    "\n",
    "    def getclipdata(self):\n",
    "        \"\"\"\n",
    "        返回一个mask处理好的矩阵\n",
    "        :return:\n",
    "        \"\"\"\n",
    "        value = self.nc_target_data.copy()\n",
    "        for i in tqdm(range(self.time_data.shape[0])):\n",
    "            temp = value[i, :, :]\n",
    "            temp[~self.mask_matrix] = np.nan\n",
    "            value[i, :, :] = temp\n",
    "\n",
    "        return value\n",
    "\n",
    "\n",
    "nc_mask = GetMask(geopandas=shp_data, nc_data=test_nc, nc_variable='temp', lat_variable='lat',\n",
    "                  lon_variable='lon', time_variable='time')\n",
    "\n",
    "nc_mask.run()"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "100%|██████████| 366/366 [00:00<00:00, 3119.49it/s]\n"
     ]
    }
   ],
   "source": [
    "clip_test = nc_mask.getclipdata()"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "markdown",
   "source": [
    "## 测试一下这个函数"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%% md\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "outputs": [],
   "source": [
    "Lon_data, Lat_data = np.meshgrid(target_lon, target_lat)\n",
    "\n",
    "plot_data = pd.DataFrame({'lon': Lon_data.flatten(),\n",
    "                          'lat': Lat_data.flatten(),\n",
    "                          'mask': nc_mask.mask_matrix.flatten()})"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "outputs": [],
   "source": [
    "% matplotlib inline\n",
    "fig, ax = plt.subplots(figsize=(7, 6), dpi=150)\n",
    "shp_data.boundary.plot(ax=ax, color='black')\n",
    "ax.grid()\n",
    "\n",
    "plot_data_in = plot_data.loc[plot_data['mask']]\n",
    "ax.scatter(plot_data_in['lon'], plot_data_in['lat'], s=0.1)\n",
    "\n",
    "plot_data_out = plot_data.loc[~plot_data['mask']]\n",
    "ax.scatter(plot_data_out['lon'], plot_data_out['lat'], s=0.1, c='red')"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Using matplotlib backend: MacOSX\n"
     ]
    },
    {
     "data": {
      "text/plain": "<matplotlib.colorbar.Colorbar at 0x7f97542418e0>"
     },
     "execution_count": 54,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "%matplotlib\n",
    "fig, ax = plt.subplots(figsize=(7, 6), dpi=150)\n",
    "\n",
    "im = ax.imshow(clip_test[1, :, :][::-1, :], cmap=plt.cm.get_cmap('RdYlBu'))\n",
    "\n",
    "fig.colorbar(im, orientation='vertical')"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 55,
   "outputs": [],
   "source": [
    "with nc.Dataset('./结果/test11152.nc', mode='w', format='NETCDF4_CLASSIC') as ncfile:\n",
    "    # 创建维度\n",
    "    lat_dim = ncfile.createDimension('lat', 360)  # latitude axis\n",
    "    lon_dim = ncfile.createDimension('lon', 720)  # longitude axis\n",
    "    time_dim = ncfile.createDimension('time', None)\n",
    "\n",
    "    # 创建变量\n",
    "    lat = ncfile.createVariable('lat', np.float32, ('lat',))\n",
    "    lat.units = 'degrees_north'\n",
    "    lat.long_name = 'latitude'\n",
    "\n",
    "    lon = ncfile.createVariable('lon', np.float32, ('lon',))\n",
    "    lon.units = 'degrees_east'\n",
    "    lon.long_name = 'longitude'\n",
    "\n",
    "    time = ncfile.createVariable('time', np.float64, ('time',))\n",
    "    time.units = 'days since 1988-01-01 00:00:00'\n",
    "    time.long_name = 'time'\n",
    "\n",
    "    temp = ncfile.createVariable('temp', np.float64, ('time', 'lat', 'lon'))  # note: unlimited dimension is leftmost\n",
    "    # temp.units = 'K' # degrees Kelvin\n",
    "    temp.standard_name = 'air_temperature'\n",
    "\n",
    "    # 写入变量\n",
    "    lat[:] = target_lat\n",
    "    lon[:] = target_lon\n",
    "    time[:] = raw_time\n",
    "    temp[:, :, :] = clip_test"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "outputs": [],
   "source": [],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}